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We consider the dynamics of a one-dimensional system consisting of dissimilar hardcore interact- 
ing (bouncy) random walkers. The walkers' (diffusing particles') friction constants £ n , where n labels 
different bouncy walkers, are drawn from a distribution g(£, n )- We provide an approximate analytic 
solution to this recent single-file problem by combining harmonization and effective medium tech- 
niques. Two classes of systems are identified: when g(£ n ) is heavy-tailed, g(£n) — £n ~ a (0 < a < 1 ) 
for large £ n , we identify a new universality class in which density relaxations, characterized by the 
dynamic structure factor S(Q,t), follows a Mittag-Leffler relaxation, and the the mean square dis- 
placement of a tracer particle (MSD) grows as t 5 with time t, where 8 = a/(l + a). If instead g is 
light-tailed such that the mean friction constant exist, S(Q, t) decays exponentially and the MSD 
scales as t 1 ^ 2 . We also derive tracer particle force response relations. All results are corroborated 
by simulations and explained in a simplified model. 



I. INTRODUCTION 

The staggering motion of a lone drunk has, since the 
classic correspondence [l], Q between Karl Pearson and 
Lord Rayleigh in 1905, been a well established metaphor 
for a random walk. This view was embraced in a cel- 
ebrated work Q by Michael E. Fisher in which he dis- 
cussed statistical aspects of many drunks. In particular, 
Fisher introduced the notion of harmless drunks for iden- 
tical non-interacting walkers, vicious drunks representing 
the motion of antagonists, and bouncy walkers exhibiting 
hardcore repulsive interactions in one dimension. 

The bouncy walker problem, which we here revisit, 
finds a number of experimental realizations: transport 
in microporous materials [H-H (e.g. zeolites), colloidal 
systems |7[ , molecular sieves [8( and biological pores [9j . 
Hard core repulsion of binding proteins diffusing along 
DNA was recently suggested to be important in tran- 
scription. [lOj 

On the theoretical side, the many-body problem of 
identical bouncy walkers in one dimension, also referred 
to as single-file diffusion ljj or symmetric exclusion [l2j , 
has received considerable attention since it was intro- 



duced by Harris 13] in 1965. The main result for thermal 
initial conditions is that a tracer particle exhibits sub- 
diffusion [n], 

EMU; the tracer mean square displace- 
ment (MSD) is proportional to t 1 ^ 2 even though the col- 
lective behavior is diffusive (the structure factor S(Q,t) 



decays exponentially jl7j])- In this study, we ask: what 
are corresponding expressions for the MSD and S(Q,i) 
for a system of dissimilar bouncy walkers? We show that 
this question can be answered in closed form, using har- 
monization [l8j and effective medium approaches [20l - 
l22j ] . The analytic predictions are corroborated by exten- 
sive simulations. Our results extend, in ways described 
throughout the study, a very limited number of previous 
studies [23-26] of bouncy walkers of different type. 



II. EQUATIONS OF MOTION 

We consider strongly overdamped motion of Brownian 
particles, in an infinite one dimensional system, interact- 
ing via a two-body short-range repulsive potential. The 
nearest neighbor potential V(\x n (t)— x n < (i)|), where x n (t) 
is the position of the nth particle, has a hard-core part 
which excludes particles from overtaking each other. The 
Langevin equations of motion are thus 



x n ,(t)]+r) n (t) + f n {t), (1) 
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where a dot denotes time derivative, f = —dV/dx n , r/ n (t) 
is a Gaussian zero- mean noise, <t] n (t)> — 0, with corre- 
lations that are determined by the fluctuation-dissipation 
theorem [13] to be <r] n (t)r) n > (t')> = 2k B T( n 6{t-t')5 n ^ , 
where ks is the Boltzmann constant and T the temper- 
ature. f n {t) is an external force and particle n has the 
friction constant £„. In Sec. IIVI where the system's col- 
lective behavior is considered we take f n (t) = 0. At the 
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Harmonization 



FIG. 1: Two mappings bring the bouncy walker problem (top) 
on a tractable form. Harmonization maps the system onto 
that of harmonically coupled, Rouse chain type beads 
(middle). The corresponding Eq. ((2j has the same form as 
the equation appearing in resistor-capacitor network theory 
where x n represent the node potentials, £ n are capacitances, 
k is the conductance of the resistors and 7/ n + /n is an external 
current entering the circuit at node n (bottom). This analogy 
allow us to derive an effective medium Eq. Q of the type 
originally used to calculate conductances of binary metallic 
mixtures [2(|] and later extended to resistor networks [2l[ and 
applied to diffusion in random media [2^ |. 



end of Sec. IVl and in Sec. [VI] we assume f n (t) to be 
a static force turned on at t = or an oscillating force 
acting only on the tracer particle. 

The dissimilarity of the random walkers enter through 
their different friction constants, £„, which here are as- 
sumed to be identically distributed random variables 
taken from a probability density £>(£«). Eq. (JTJ) , the 
noise correlation and g(£, n ) defines the dissimilar bouncy 
walkers problem which we deal with by using two non- 
equilibrium statistical physics methods outlined below 
and elaborated in appendix 1X1 

We complement all analytic results by stochastic sim- 
ulations for hardcore interacting particles in a box with 
reflecting boundary conditions, detailed in appendix [G] 
The box size is chosen sufficiently large such that bound- 
ary effects are negligible for the center particle within the 
duration of the simulation run. 

A prototypical system where Eq. is applicable is 
protein diffusion on DNA molecules [lOj. where the het- 
erogeneity in friction constants originate from the dif- 
ferent proteins varying binding strength to the DNA 
molecule. 



III. HARMONIZATION AND EFFECTIVE 
MEDIUM THEORY 

Our first step towards bringing the dissimilar walker 
problem on a tractable form is to generalize the harmo- 
nization approach (l8j to the case when the particles have 
different friction constants. In the long time limit this al- 
lows us to map the equations of motion ([I]) onto that for 
a linear chain of interconnected springs 
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dx n (t) 
dt 



k [x n+1 (t) + - 2x n (t)]+r] n (t)+f n (t), 

(2) 

corresponding to harmonically coupled beads, see Fig. 
[T] The effective nearest neighbor spring constant k is 
calculated by demanding that the change in free en- 
ergy for small particle displacements in the original sys- 
tem equals the free energy change in the spring sys- 
tem. For hardcore interacting particles of size b (used in 
our simulations) this harmonization procedure yields [18[ 
k = p 2 kBT(l — pb)~ 2 , where p is the particle density. The 
harmonization approach relies on the assumption that lo- 
cal equilibration is faster than tracer particle dynamics, 
for sufficiently long times. We find self-consistcntly that 
this holds since, as we will see, the MSD for a single par- 
ticle is proportional to t s with S < 1/2 [see Eq. (fTTj) ]. i.e., 
the particle cross a distance of length L in a time on the 
order of L 2 / s . At long distances this is indeed slower than 
the corresponding density relaxation time which scales as 
L 1 / 6 [see Eq. ©]. 

Second, applying an effective medium approach [20l — 
l22j | to the harmonized equations yields a set of gener- 
alized Langevin equations containing a time dependent, 
but n independent, memory kernel £ e ff(i)- This method 
is based on an analogy of our system with that of resis- 
tor networks, see Fig. [TJ In the continuum limit (in the 
long-wavelength limit n can be treated as a continuum 
variable [18|) the effective medium harmonized equations 



Z eS (t-t')x n (t')dt' 



"^^-Hjf (*)+/«(*)• (3) 



The lower integration limit being — oo physically corre- 
sponds to the system dynamics being "turned on" in 
the infinite past, thereby effectively bring the system 
to equilibrium at any finite time t. As shown below 
£,cff(t — t') = for t' > t, as required by causality. 
The fluctuation-dissipation theorem enforces a relation 
between the memory kernel and the effective noise rf^ (t) 

«(t)ftf>tf)) - k B TCM\t - t'\)S(n - n'), (4) 

where the brackets (..) represent an implicit average over 
quenched friction constants, besides averages over dif- 
ferent realization of the noise and random initial posi- 
tions. We label this average the heterogeneity-averaged 
case [23j j . which is contrasted by the non-averaged case 



3 



(represented by <..>) where an average over the proba- 
bility density of friction constants is not performed. In 
the simulations for the non-averaged case the same £ n 's 
are used when averaging over thermal noise (i.e., for each 
simulation run). For the heterogeneity averaged case we 
draw new friction constants whenever we make a new 
initial particle positioning. 

The memory kernel, £ ff(f)j in Eq. ([3]) is determined 
by imposing a self-consistency criterion (see appendix^} ■ 
We identify two classes of systems: For light-tailed (LT) 
systems with a probability density g(£) such that the 
mean friction £ = J £p(£)d£ exists, we obtain, in the 
long time limit, a memoryless kernel 



(5) 



where S(z) is the Dirac delta-function. Thus, LT sys- 
tems are in the same universality class as that of identi- 
cal bouncy walkers |l8]. In contrast, we find that heavy- 
tailed (HT) systems, where g(£) — A^ 1 ^ for large £, 
with < a < 1 (A is a normalization constant), belong 
to a new universality class where the memory kernel has 
a power-law decay with time 



-2(5 



£eff(*) ~ X 



r(i - 25) 



0(f) 



(6) 



where T(z) is the gamma function, 9(z) is the Heaviside 
step function, X = (^n) 1 - 26 (Att / sm{5n / (1 - S)]) 2 ^ 
and 



6 = 



1 + a 



(7) 



We point out that for HT systems it is only the tail of 
the distribution [the detailed structure of q(£) for small 
£ is not important] which determines the long-time dy- 
namics. The equations ©-© in this section allow us to 
calculate explicit long-time expressions for observables in 
dissimilar bouncy walker systems [i.e., systems described 
by Eq. © . 



relaxations are thus exponential as for a system of harm- 
less drunks 28] or identical bouncy walkers [25[ . For HT 
systems we find a collective behavior which is drastically 
different: 



S(Q,t) ~ S(Q,0)E 2S (-\ c Q't 



2,2S\ 



Q^O (9) 



where E a (z) is the Mittag-Leffier function [29[ and A c = 
K /(p 2 x) is a generalized collective diffusion constant. 
Here, the density relaxations exhibit anomalously slow 
power-law decay with time. The anomalous decay of 
S(Q,t) is illustrated in Fig. [^together with simulation 
results for the non-averaged case. 

From Eq. © and Onsager's regression hypothe- 
sis [13] it follows that perturbations of the concentra- 
tion c(X,t) = (%2 m 6(X — x m (t))) around its equi- 
librium value p decay according to the fractional dif- 
fusion equation dc(X,t)/dt = \ c D l t ~ 2S d 2 c{X , t)/dX 2 
where qD\ ~ 2S is the fractional Riemann-Liouville opera- 
tor. This equation describes also the subdiffusive mo- 
tion of non-interacting continuous time random walk- 
ers (CTRW) with a power-law waiting time density [29l ]. 
However, the nature of the process is here very different. 
In particular, our system has a stationary fluctuating 
equilibrium state, since the underlying many-body dy- 
namics is Markovian with a well-defined equilibrium (two 
point correlation functions only depend on the time dif- 
ference, i.e., the system do not age). In contrast, CTRWs 
do not have a stationary state (their two point correlation 
functions age). Therefore, to the best of our knowledge, 
dissimilar bouncy walkers are the first example of a sta- 
tionary physical system obeying (to a good approxima- 
tion) the fractional diffusion equation for density relax- 
ations. The fact that our HT bouncy walkers belong to 
a different universality class than CTRW systems is also 
captured by the difference in tracer particle behavior. 



V. TRACER PARTICLE DYNAMICS, 
HETEROGENEITY-AVERAGED CASE 



IV. COLLECTIVE BEHAVIOR 

Let us now consider the collective behavior of dissim- 
ilar bouncy walker systems, using the analytic approach 
in Sec. IIIII and the simulation scheme from appendix 
[Gl A quantity capturing the collective motion and easily 
accessible in experiments [28j is the dynamic structure 
factor S(Q,t). For translationally invariant systems we 
use S(Q,t) = J2 n {e iQ ^ Xn ^- xo( - ^) where the summation 
index runs over all particles. For LT systems we find from 
Eq. (see appendix [C]) , in the limit of small wavevec- 
tors Q ^ and long times, that 

S(Q, t) ~ S(Q, 0) exp{-D c Q 2 t) (8) 

with the static structure factor S(Q, 0) = ksT p 2 j k and 
the collective diffusion constant D c = n/(p 2 ^). Density 



We now address tracer particle dynamics in the ab- 
sence of an external force (see appendix [D] for details) 
using Eqs. (El)-© . For LT systems we find the MSD 



(6x 2 T (t))~k B TJ—-cxt 1 / 2 , (10) 

for particle T(= 0). Thus, the MSD for LT systems 
take the same form as for identical walkers [18J, but with 
the friction constant for the homogeneous case replaced 
by the mean friction constant The same f 1 / 2 -scaling 
was recently obtained for dissimilar walkers on a lattice. 
[H HI] A very different MSD exponent is found for HT 
systems where instead 

(5x 2 T (t))~^-f—oct s (11) 
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FIG. 2: Dynamic structure factor. Collapse-plot for S(Q,t) 
for different wavevectors Q as a function of time, illustrating 
anomalous decay. The black solid curve is the analytic pre- 
diction, Eq. (0. Simulations are non-averaged results from 
32000 noise realizations with N = 501 particles in a system of 
length L=2501b (6 is the particle size) with a = 0.5. The av- 
erage diffusion constant is denoted D. Note that, in contrast 
to tracer particle observables (see inset of Fig. [3} S(Q,t) is a 
sum over all particles and therefore display self-averaging. 

with 5 < 1/2 [see Eq. 0] indicating ultra-slow dynamics 
for the tracer particle. Heterogeneity-averaged simula- 
tions show excellent agreement with this result, see Fig. 
[3j The MSD exponent S — a/(l + a) has also recently 
been derived for lattice systems in a rather technical 
mathematical study by Jara [26| and obtained through 
scaling arguments in Ref. l2~il Besides providing a simpli- 
fied derivation for continuum systems, we extend Jara's 
result by also obtaining an explicit expression for the 
MSD prefactor. The probability density function (PDF) 
for the tracer particle position is a Gaussian (since the 
noise is taken from a multi-variate normal distribution) 
with a width described by Eq. (JTTJ) , see appendix [G"l for 
corroborating simulations. 

In the presence of a external force on the tracer 
particle we calculate force-response relations (see Ap- 
pendices [E] and |F|. For the case of a static force 
fn(t) — S n ,oFo@(t) with magnitude Fq, we show by ex- 
plicit calculation that the generalized Einstein relation 
(x T (t)) f = (5a%-(t))F /(k B T) is satisfied for both LT 
and HT systems. The bracket (..)/ denotes an aver- 
age in the presence of f n (t)- For an oscillatory force 
fn(t) = SnflFo cos(woi) with angular velocity ljq we get 

{x-f{t))f = [F /{2^J k^usq)] cos(woi — </>) where (j) = 7r /4 
for LT systems. For HT systems we find 

F 

(x T (t)) f = - — g cos(w * - 4>) (12) 
with the non-trivial phase-shift <f> — Stt/2 < 7r/4, 



FIG. 3: Tracer particle mean square displacement. MSD as 
a function of time for the heterogeneity-averaged and non- 
averaged (inset) cases. The straight solid lines are the asymp- 
totic analytic prediction, Eq. (jTTj) , and the dashed lines 
are the result for an independent random walker. In the 
heterogeneity-averaged case the simulations are averages over 
2400 realizations, and the center particle is taken as tracer 
particle. We used a box size L — 100016 (b is the particle 
size) and N — 1001 particles for a — 0.7 and a = 0.5. For 
a = 0.3 (a = 0.1) we had N = 501 and L = 5001& (N = 251 
and L — 25016). The average diffusion constant is D. Inset: 
Non-averaged MSD simulations for a = 0.5 and 800 noise 
realizations for each curve. Notice the spread around the an- 
alytic result. 

which compactly quantifies the subtle interplay between 
bounciness and the degree of friction dissimilarity. Eq. 
(1121) is compared to simulations for the heterogeneity- 
averaged case in Fig. 01 

VI. SELF-AVERAGING 

Let us comment on the difference between the non- 
averaged and heterogeneity-averaged cases. From simu- 
lations for both cases we found that S(Q, t) and the MSD 
are self-averaging quantities for the LT class, as has been 
previously shown for lattice systems [23j, |25( . For HT sys- 
tems S(Q, t) self-averages while the MSD does not (see 
Figs. [5]and[2]). In order to understand this in more detail, 
we put forward a simplified model for tracer particle dy- 
namics in the non-averaged case: under the influence of a 
constant force, Fo, we write £ s totd<8xi-> p / dt = Fq. We 
use brackets <..> to denote averages over thermal noise 
and random initial positions for the non-averaged case. 
The quantity £tot represents the total friction experienced 
by the particle which is not only its own, but also that 
of an approximate number TV of its neighbors on which 
it is "pushing" , i.e. £ tot = Co +£i + 6 + ■ • ■ + £,n ■ Within 
the harmonization approach [18| we have approximately 
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FIG. 4: Mean displacement for a tracer particle in the 
presence of a time-varying force. Simulations for the aver- 
age displacement of a tracer particle as a function of time 
within a single period of a harmonically oscillating external 
force, heterogeneity-averaged case. Parameters are a — 0.5, 
TV = 501, L = 50016 (b is the particle size), a = 0.002 and 
ujo/(2tt) = 10 _J . The average is over last 100 periods of the 
force for 1280 noise realizations each lasting 200 periods. The 
blue solid curve is the analytic prediction, Eq. (|12[) . Notice 
the excellent agreement between theory and simulations. 



that the effective spring constant for TV serially connected 
springs is ft/TV, i.e., we have F — (k/N)<Sxj->f , from 
which TV = k<Sxj->f /Fq; the larger the displacement 
the larger is the number of particles contributing to £ t ot ■ 
Using this result and multiplying and dividing the equa- 
tion of motion above by TV 1 /" we arrive at 



6+6 



N 



TVV« 
C 



(k<Sx t >f q /F ) 



1/a d<Sx r > Fn 



dt 



where £ is a random variable with an TV-independent 
distribution V(C) f° r N ~* 00 • Solving the equation 
above and employing the fluctuation-dissipation theo- 
rem in the form of the generalized Einstein relation, 
k B T<5xr>F /F = <Sx^->, we find the MSD for the 
non-averaged case: 



<5x 2 r> oc 



-t" 



(13) 



In agreement with the non-averaged case simulations (in- 
set Fig. 13]), the MSD prefactor is a random variable £. 
Thus, HT systems do not show self-averaging, rather each 
new realization of £„'s gives a different fluctuating pref- 
actor for the MSD (but with the same scaling with time) . 
The very simplistic argument here predicts that for HT 
systems tp(Q is a one-sided Levy stable distribution [2Tj| |. 
but further work is needed to pinpoint the correct func- 
tional form for ■ip(C). We point out that the non-averaged 



HT case has not been investigated in previous studies 
Using the fact that ((~ s ) oc A^ 1 " 5 ) for the power- 
law distribution of friction constants considered here, our 
simplified model recovers Eq. (|11[) , up to a dimension- 
less a-dependent prefactor for the heterogeneity-averaged 
case. For LT systems (a = 1) the central limit theorem 
gives tp(() — <5(C — as TV — !> oo, i.e., the argument above 
predicts that LT systems possess the self-averaging prop- 
erty as indeed found in our simulations, pfjj 



VII. SUMMARY AND OUTLOOK. 

We investigated the dynamics of a one-dimensional sys- 
tem consisting of dissimilar bouncy random walkers, with 
friction constants drawn from a probability density. Two 
classes of systems were identified: those with heavy-tailed 
(HT) probability densities in which density relaxations, 
characterized by the dynamic structure factor S(Q,t), 
follows a Mittag-Leffler relaxation. The mean square 
displacement of a tracer particle (MSD) grows as t 5 with 
time t, where 5 < 1/2, and we find a phase-shift <p = Sir/ 2 
in the force-response relation for such systems. For light- 
tailed (LT) probability densities, S(Q,t) decay exponen- 
tially, the MSD scales as t 1 / 2 and a phase-shift <j> = it/ A 
in the force response is obtained. We also introduced a 
simplified model which allowed us to address the prob- 
lem of self- averaging for LT and HT systems. All results 
were corroborated by extensive simulations. 

Prior to this study, dissimilar hardcore interacting dif- 
fusing particles have been investigated only in a few pub- 
lications, see Refs. 12314261 Using harmonization and ef- 
fective medium approaches we here extended previous 
results by providing: (i) an explicit expression for the 
structure factor S(Q,t); (ii) an explicit prefactor for the 
MSD for HT systems; (iii) force-response relations for 
a tracer particle; (iv) insights concerning the difference 
between non-averaged and heterogeneity-averaged cases 
for HT systems; (v) extensive simulation results. We 
further provided a much simplified framework, Eqs. ©- 
©, within which further observables, such as correlation 
functions between particles (see Ref. [l8|), can be calcu- 
lated. 

Solvable many-body problems have served as impor- 
tant models for real systems over the last century. The 
problem of dissimilar bouncy walkers is expected to find 
applications in systems such as passive diffusion of pro- 
teins along DNA or particles diffusing in nanochannels. 
Further, often the dynamics of macromolecules in living 
cells is found to be subdiffusive [3l|. The techniques de- 
veloped in this study may be of use also in such higher 
dimensional heterogeneous systems. 
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Appendix A: Effective medium approach. 

Our second step in deriving a set of manageable equa- 
tions for the dissimilar bouncy walker problem is to in- 
troduce an effective medium approximation known for in- 
stance from the theory of resistor networks [32| , see Fig. 
[TJ In the language of single-file motion the approxima- 
tion consists of replacing the different frictions in Eq. 
(12) with an effective friction, or memory kernel, £ c ff(i), 
which is identical for all the particles, but instead has a 
memory. The memory kernel £ e ff (t) is chosen such that if 
one takes a random particle from the harmonized system 
with friction constant £ n and place it among particles 
characterized by £ e ff(*), then the mobility of this parti- 
cle is, on average [averaged using £>(£«)], required to be 
identical to the mobility of a tracer particle in the sys- 
tem where all particles had friction £ e g(t). A detailed 
calculation of £ e S (t) is given below. 

Let us, as our starting point, consider the harmonized 
equation which we write: 

Cn^^ = « [Vn+l(t) + Vn-l{t) - 2 Vn (t)] + T) n (t) + /„(*) 

(Al) 

where we introduced y n (t) — x n (t) — np , which elimi- 
nates the average values so that <y n (t)> = 0. Our goal 
here is to find the friction kernel £ c ff(*) entering the ef- 
fective medium harmonization equation of motion 



&ff (* - t')-^^-dt' = k [y n +i(t) + y n -i(t) - 2y n (t)} 



+vT(t) + fn(t) 



(A2) 



which in the continuum limit with respect to n becomes 
Eq. ([3]) in the main text. We introduce the frequency- 
dependent generalized friction constants 



7^H 



F, 



nil— >n 



<v n (u>)> 



(A3) 



where v n (uS) = —iu)y n (uS) is the Fourier transform of the 
velocity, and F n ^ n+ \ — —n(y n+ i — y n ) is the force from 
particle n on particle n + 1 within the harmonization ap- 
proach. The quantity 7« (7n ) can be thought of as 
a generalized friction of particle n due to the part of the 
harmonic chain lying to right (left) of particle n. We here 
define the Fourier-transform with respect to time of a 
function A(t) according to: A(u>) 

inverse transform A(t) = J 



joo^ e i^A(t)dt, with 



•~ tu]t A(u})du;/(2n). Us- 
ing the generalized frictions introduced above and setting 



fn(t) = fo(t)3n,o one readily obtains, from the Fourier- 
transform of Eq. @, the mobility of particle 



_ («o(w)> 



1 



/oH £o+7^ +) H+7rM 
as well as the recursion relations 



(A4) 



7i +) H- 



1 — iu 



n>0 (A5) 



and 



In 



Cn-l+7^_l(w) 



Cn-l+7i-lH /« 



n < (A6) 



So far everything is only a reformulation of Eq. (IA1|) 
for the case f n (t) = fo(t)S n fl- We point out that we 
could equally well have reformulated Eq. (1A1[) by using 
Laplace-transformed quantities A(s) = e~ st A(t)dt; 
the corresponding recursion relations would then take the 
same form as in Eqs. (|A5[) and (IA6|) . with s = —ioj. How- 
ever, for the purpose of avoiding to explicitly compute 
averages over initial positions using Fourier-transforms 
is more convenient, and we therefore employ Fourier- 
transform techniques throughout this study. 

Let us now turn to the problem of computing £ c ff (i) hi 
the effective medium Eq. (| A2|) . We invoke the following 
procedure: (A) Fourier-transform Eq. (|A1|) . and replace 
all friction constants for n 7^ by an n independent but 
frequency dependent effective friction, i.e. we let 



Cn (1 - #n,o)feff(w) + £o^n,0- 



(A7) 



From these equations the mobility Mo(£o) for particle 
is determined. (B) We then obtain £ e g(uj) by impos- 
ing the self-consistency requirement that /io(£o) averaged 
over different realizations of £o 5 



/'o 



Mo(£o)£>(£o)d£o, 



(A8) 



is equal to the tracer particle mobility, (j, e s(u>), obtained 
from the effective medium Eq. (|A2[) . i.e., when also £q is 
replaced by £ c s(lu). 

Let us now carry out the scheme above. (A) From Eqs. 
(lA5l) . JM| and (TATl) 



7eff 



£eff + 7cff 



f 



1 [£eff +7off] /«' 



where we left arguments implicit, and 7 e g = 7n 
(independent of n) . Thus, 



(+) 



7off 



£cff 
2 



;)'( 



(A9) 



In 



(A10) 



where we used the fact that the real part of 7 e s needs 
to be positive (so the friction will dissipate energy) to 
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choose the correct root of the quadratic Eq. 
small frequencies (long times) we have 



7eff 



S). For 



(All) 



(B) From the self-consistency criterion we now deter- 
mine £ e ff- Setting n = in Eqs. (IA5I) and (|A6I) we 
obtain the generalized friction 7^ = (£ c ff + 7 G fr)/(l — 
[CeS + 7eff] /k) for particle 0. Combining this result 
with Eqs. (|A4I) . (|A8|) and the self-consistency criterion, 
A) = ^cff(w), we get: 



1 



g(£o)rf£o 

£0 + 27eff £eff + 27cff 



(A12) 



where the right-hand side is the mobility obtained 
through the effective medium equations [simply replace 

£0 -> Coff and 7^ ±} -> 7 cff in Eq. flU}]. Eq. (jAl~2)) can 
be rewritten: 



elf 



g(£o)rf£o 
Co + 27 cff 
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g(£o)d£o 
Co + 27 cff 



(A13) 



where the normalization condition 1 = Jp 00 g(£o)^£o was 
used. To summarize briefly, the effective friction £ c ff (<-"-") 
is uniquely determined by Eqs. (IA10[) [Eq. (jAllI) for 
small frequencies)] and (|A13[) for a given choice of g(£o)- 
Now we have to distinguish two cases: (LT) the distri- 
bution £?(£o) nas a finite first moment, £ = J £og(£o)d£o < 
00. In this case the £o's in the denominators on both sides 
of Eq. (|A13p can be discarded for large oj since we expect 
(and find self-consistently) that 7 e ff(w) -4 00 as w -> 0. 
For LT systems we therefore have: 



£eff(w) ~ f. 



(A14) 



which in the time-domain agrees with the result in the 
main text. A more detailed analysis shows the result 
above applies for small frequencies in the sense that 
w£/k <C 1. (HT) in the case of g(£o) not possess- 
ing a finite first moment we assume a power-law tail 
g(£o) ~ ^0 f° r large £0 (0 < a < 1). In this case the 
£0 in the denominator on the left hand side of Eq. (|A13|) 
can again be neglected at small uj. In contrast, on the 
right hand side of Eq. (| A13|) the tail of g(£o) dominates 
the integral at small w; this integral therefore needs to 
be evaluated. We then find, using Eq. (IA11[) . that Eq. 
(IA13j) becomes: 



\ —lUJ J 



-1/2 



Att /4< cff \- Q / 2 



sin(7ra) V — iu> 



V — iui ) 



(A15) 



so that for HT systems the effective friction, for small 
frequencies, decays according to a power-law: 



feff(w) 



X 



(-iw) 



1-25 



(A16) 



with a prefactor X = (4k) 1 - 26 (An / sin[57r/(l - 5)}) 2{1 - S) 
and where <5 = a/(l + a). The result above applies for 



the case that uj\^ e s(uj)\/ k <C 1. Fourier-transforming to 
the time-domain we obtain Eq. © in the main text. 

Let us finally, for completeness, also give the expres- 
sions for the tracer particle mobility /io(w). Combining 
Eqs. (|A4[) and (|A11[) and taking the u) — > limit, we get 
/x (w) ~ (-iw) 1 / 2 /[2(<) 1 / 2 ] for LT probability densities, 
and 



(-iu) 



1-8 



2(«X) 1 / 2 



(A17) 



for HT systems. 



Appendix B: Solution of the effective medium 
equations. 

Taking the Fourier-transform with respect to n and t, 
the solution to Eq. ((3]) in the main text is 



7 ? ctt (g,w) + /( g ,w) 
nq 2 - iw£ e ff(w) 



(Bl) 



where y n {t) — x n {t) — np , and the Fourier-transforms 
of with respect to n of a function B(n) is de- 



fined: B(q) = j_ 



-tqn 



B(n)dn with inverse B(n) 



f_ oo e lqn B(q)dq/(2ir). Also Fourier-transforming Eq. (@| 
of the main text leads to 

(^fawto^a/)) = {2n) 2 k B T[^)+^ s (uj')] 

x5{u + u')5(q + q'), (B2) 



which together with 



(rf tt (g, W ))=0 



(B3) 



provides a general solution to the effective medium har- 
monization equations. Using Eq. (IB1[) we find 

, _ r 00 dq r°° d£ j 00 dw r°° do/ 



e i(gn+g'n')e i ( w *+ w ' t ')^ efr ( g) a;) 7? eff (g / ) a; / )) 

[k^ 2 - iu;£ off (u;)][K<7 /2 - iu)'^ eB (u)')] 



(B4) 



so that with the help of Eq. (IB2|) and some simple alge- 
braic manipulations we obtain 



(y n (t)y n '(t')) = 



1 



^9 iq(n-n') 

2tt 

1 



2^ 



iw q 2 — iw£ c ff (w) / K 



-iui(t-t') 



(B5) 



where c.c. denoted the complex conjugate. The relation 
above will be used in subsequent sections. 



Appendix C: Dynamic structure factor 

The dynamic structure factor is defined: S(Q,t) — 
Enm( e '*°' f ' _Im ' ll ")/^i where the sums extend over 
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all particle labels [28j. For an infinite, translationally 
invariant system one can immediately perform one of the 
sums and get 



was bounded from above we have then to first order in 



S(Q,t)= J2 



o iQ[x n (t)-x (0)]\ 



(CI) 



Furthermore, assuming that the process is Gaussian one 
knows the characteristic function and can therefore eval- 
uate the average over noise to find 



S(Q,t)= e iQ ^- Q2 ^ 



(t)/2 



(C2) 
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where we have introduced the average and variance 

M*) = ( x n{t) - x {0)) = np-\ (C3) 
= (M*) - *o(0) - M„(<)] 2 ) 
= ([y n (t)-y (0)} 2 ). (C4) 



The variance can be evaluated using Eq. (JB5J) as 

iqn—iujt \ 

2ir ^ ' 
AbT 1 



+ c.c. (C5) 



\ Kg 2 K(? 2 /^ ff (w) — 



At this point we introduce a collective diffusion con- 
stant A c by writing £ c ff (<-"-") = K A f T 1 (— itu) 25 ^ 1 / p 2 . For 
LT systems we see from Eq. (|A14[) that 5 = 1/2 and 
A c = k/(p 2 £), and for HT systems Eq. (|A16|) gives 
A c = k/(p 2 x)- Now we perform the w-integration by 
changing variable — iu) — > s and using that a definition of 
the Mittag-Leffler function E a (-) is 



E a {-et a )e- st dt 
o s a +9 



(C6) 



to find 



->,,•, i / dq k B T , 2\ „2ui2(5\ 



27T K(/ 2 

This result can also be written 
k H T 



(l-e^ n E 2S (-p'X c q'\tn). (C7) 



^(*) = ^- (N + Px/A:|t|^(n/(pV^:|t| 5 ))) , (C8) 
where we have introduced the function 

h(z) = 2 /" g cos(«z)[l - ^(-g 2 )]/^. (C9) 

Note that is bounded from above by /i(0). 

To get further, we are now going to take the hydrody- 
namic limit of Q — > and t — > oo. In taking this limit 
we keep Qt^ fixed, where /3 is an exponent which will be 
chosen such that we get a non-trivial result. If we make 
the ansatz that > 6/2 then Q 2 t s —> and since h(z) 



(CIO) 

X U - \Q 2 ^pVMt\ 5 h(n/{p^\ c \t\ & )^ . (Cll) 

If we insert this into the expression for S(Q,t) and take 
the continuum limit by replacing the sum with an 
integral J dn then the expression for S(Q, t) can now be 
evaluated. For (3 < 6 one gets S(Q, t) ~ while for (3 > 5: 
S(Q,t) ~ ksTp 2 /K. The non-trivial result appears for 
(3 = 8 where one finds 



S(Q,t) 



k B Tp 2 



E 2 s(-KQ' 2 \t\ 25 ) 



(C12) 



This is the expression for S(Q,t) used in the main text. 

As mentioned in the text, the Mittag-Leffler decay 
of the structure factor implies, by Onsager's regression 
hypothesis [301 ] . that perturbations of the concentration 
c(X, t) = (J2 n S(X — x n (t)) decay according to the frac- 
tional diffusion equation (FDE). To see this note that 
the structure factor is the Fourier transform of density- 
density correlations. Onsager's regression hypothesis (al- 
ternatively the fluctuation-dissipation theorem) implies 
that perturbations in the density decay in the same way 
as correlations. Thus we have for the Fourier transform 
c(Q,t) = J dX e~ l ® x c(X, t) that it decays according to 
c(Q,t) = c(Q,t = 0) x E 2S {~\ c Q 2 \t\ 2S ) (Q > 0). This 
implies that it obeys a fractional relaxation equation 

J^c(Q, t) = -X^oDl^ciQ, t), (C13) 

where the fractional Riemann-Liouville operator qD" is 
defined for < a < 1 by 



D?c(Q,t) = 



1 



d 



r(i 



rt dt , c(Q,t') 
a) dt J (t - try 



(C14) 



Taking the inverse Fourier transform of Eq. (IC13I) one 
arrives at the FDE 



8 



„ d 2 

r"'"' _ " cu "' dX 2 



f t c(X,t) = X c0 D 1 -" i 
given in the main text. 



c(X,t). 



(C15) 



Appendix D: Tracer particle mean square 
displacement. 

Let us now derive an expression for the MSD in the 
absence of an external force. We have 

(Sx r (t) 2 ) = ([ yo (t) - y (0)] 2 ) 

= 2[(y (0) 2 ) - (y o (t)y o (0))} (Dl) 

due to stationarity. Here T = has been chosen with- 
out lack of generality. Using Eq. (|B5[) and the integral 
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FIG. 5: Collapse-plot for the tagged particle probability den- 
sity function at different times. The solid black curve is a zero- 
mean Gaussian with a width given by Eq. (JTTJ . The symbols 
correspond to heterogeneity-averaged simulations over 8000 
noise realizations with N = 501 particles in a system of length 
L=5001b (b is the particle size) with a = 0.7. 



J^° oc (q 2 + C) 1 dq = n/VC for Re[C] > 0, we can express 
the time-derivative of the MSD according to: 



— (ox 
dt x 



k B T 
2^ 
1 



due 



(D2) 



Finally, rewriting the expression above as a sine- 
transform, inserting Eqs. (|A14I) and (|A16[) and integrat- 
ing [using the fact that (Sxt(0) 2 ) — 0] we recover the 
MSDs given in the main text. 

The noise is assumed to be Gaussian distributed, and 
therefore the tracer particle PDF is a Gaussian with zero 
mean and a width determined by the MSDs given in Sec. 
fVl Fig. [5] compares heterogeneity-averaged simulation 
results for the PDF to this prediction, finding excellent 
agreement. 



Appendix E: Tracer particle dynamics in the 
presence of a time-varying external force. 

To test the effective medium harmonization approach 
further we have considered the response of a tagged 
particle to an external force. For this purpose let us 
assume that a harmonically oscillating force f n (t) — 
S n fiFo cos(woi) acts on particle 0, where ujq is the angu- 
lar velocity and Fq the force amplitude. From the formal 
solution Eq. (jBlJ and Eq. (JB3) we find: 



y{q,u) 



nF [S(uj - cl>q) + 5(u + uj )} 
nq 2 - iuj£ cS (u}) 



(El) 



Further, by performing the inverse Fourier-transform 
with respect to q, and setting n = 0, we get 



(yo(t))f = f^Re 



a iui t 



-iwo^off^o)/^] 1 / 2 



(E2) 



where Re{..} denotes the real part and the subscript in- 
dicates an average in the presence of a force. For a LT 
[HT] system we use the friction constant in Eq. (|A14p 
[Eq. (TAT6)) ] in the expression above. Switching to the real 
tagged particle coordinate xj- assuming (xf(t = 0)) = 

we get (xr(t)) f = [F /(2^/ k^cjq)] cos(wo^ — 4>) where 
<\> = 7r/4 for LT systems. For HT systems we find 



Fq 

(5x T (t)) f = — s cos(w i - 0) 
with the non-trivial phase-shift <fr = Stt/2 < 7r/4. 



(E3) 



Appendix F: Generalized Einstein relation. 

Suppose a static force Fq is applied to particle at 
time t = 0, i.e. f n (t) = F 5 n fi6{t). From the fluctuation- 
dissipation theorem [23] , we have that the following 
simple relation holds between the tagged particle mean 
square displacement in the absence of the force (&c^(i)), 
and the mean displacement in the presence of the force 
(5x T (t)) F 



(dx T (t)) 3 



Fq 
2k B T 



(Sx 2 T (t)) 



(Fl) 



which is referred to as the generalized Einstein rela- 
tion. A consistency check of the effective medium- 
harmonization formalism is thus to check that this re- 
lation holds. From Eq. (|B1|) we straightforwardly com- 
pute the average velocity (in frequency space) : (v n (u))) = 
—iuj(y n (uj)) for the type of force considered here. Fur- 
ther, by Fourier-inverting and setting n = T = 0, and 
using Eq. (|A16I) we find 



(Mm 



Fq t 



6-1 



9(t) 



(F2) 



for HT systems. By integrating the expression above over 
time and comparing with the MSD in the main text we 
find that Eq. (IF1[) is satisfied for t > 0. A similar calcu- 
lation for LT systems shows that indeed the generalized 
Einstein relation is satisfied also then. 



Appendix G: Simulation scheme 

We consider N hardcore particles of linear size b in a 
box of length L (the box extends from —L/2 to L/2). 
At random times the particles performs hops with the 
jump length taken from a Gaussian distribution of width 
a. Particle n (n = 1, N) has a hop rate q n , where q n is 
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related to the friction constant by q n = 2fcsT/(^„a 2 ). 
We are interested in the stochastic dynamics of the center 
positions X n of all particles, when we have the single-file 
constraint X n+ \ — X n > b at all times. Also, particles 
cannot move out of the box: X\ > — (L — b)/2 and Xn < 
(L — b)/2. Our numerical approach detailed below, builds 
on a previous scheme for lattice simulations called the 
trial-and-error Gillespie algorithm [25j . Our continuum 
version of that algorithm is: 

1 . Generate the cumulative sums of the free total rate 
constants 

Po = 0, 

71 

Pn = 51 9m, n=l,...,N. (Gl) 

m— 1 

2. Generate a random initial configuration of particle 
positions X n (0) (n = 1,.,.,N). This step requires 
some care, see next section. 

3. Set the time equal to zero, t = 0. 

4. Draw a uniform random number r (0 < r < 1). 
Construct a waiting time r according to: 

r = — iog(-), (G2) 
Pn \rj 

i.e. t is taken from an exponential distribution 
<j>(r) oc exp(— r Yl n =i Update the time i — > 
t + r. 

5. Draw a new uniform random number r (0 < r < 1). 
A trial particle n is determined by the n which 
satisfies 

Pn-i < rp N < p n (G3) 
i.e., particle n is chosen with a probability oc q„. 

6. Draw a random attempt jump length I from 
the Gaussian distribution of width a: P(l) — 
(27m 2 )- 1 / 2 exp[-(l - fi) 2 /(2a 2 )}. For the case of a 
harmonically oscillating external force on the tracer 
particle f n (t) — S n fiFo cos(u>o£)j M = ao cos[u;oi]) 
for the tracer particle with ao = Fo/foqo), other- 
wise /i = 0. 

7. Loop over the neighboring particles in the direction 
of the attempted jump, and find the number M of 
particles which have their center positions within a 
distance I — b of their previous neighbor. 

8. Draw a uniform random number r and move par- 
ticle n and the M neighbor particles a distance / if 

r < W Emt=n £™ ( tne sum ran g es f rom n t° n+M 
if the attempted jump is to the right, otherwise the 
range is n — M to n) , and there is no wall prevent- 
ing the move. If the move occurs then update the 
particle positions. 



9. Return to step (4). 

The scheme above produces a stochastic time series 
X n (t), and if steps (2)- (9) are repeated many times, en- 
semble averages can be computed. Note that step (1) 
does not have to be repeated for each simulation run 
(each realization of the noise) for the non- averaged case, 
i.e. if the free total hop rates are the same for all simula- 
tion runs. For the heterogeneity-averaged case, where we 
assign new hop rates for the particles for every new simu- 
lation run, step 1 has to be repeated for each one of them. 
An efficient method for performing the search for the n 
satisfying the inequality in Eq. (|G3p was presented previ- 
ously (25(. We here point out that in our previous study 
[25J we instead of Eq. (|G3[) used p n -i < rpiy < Pn- This 
criterion for choosing a trial particle causes the algorithm 
to get stuck for the (exceptional) case that p n -\ = rpN 
if combined with the search algorithm in appendix B of 
Ref. [25|. We are grateful to Karl Fogelmark for pointing 
this out to us. 

In step 2 in our algorithm, we set the initial positions 
randomly (with the possible additional constraint of a 
fixed position for a tracer particle). For point-particles 
this present no problem - for completely random initial 
positioning one simply needs to draw N uniform random 
numbers between —L/2 and L/2, which are then sorted 
and assigned as initial positions. For finite size parti- 
cles one must make sure that there is no overlap of the 
particles. The same procedure as used above, with the 
additional step to discard configurations where two or 
more particle overlap, is very inefficient for large iV. We 
instead use the following mapping between the finite-size 
particle problem [initial positions X n (0), n = l,...,N] 
and a point-particle problem [initial positions X n (0)] 

nam 

L = l + Nb, 
X n (0) = X n (0)+nb~^p-b. (G4) 

Step 2 in the simulation scheme then becomes: draw ran- 
dom numbers between —1/2 and 1/2 which, after sorting, 
gives X n (Q). Then use the relation above to determine 
X n (0). This procedure is used in the simulations provid- 
ing the dynamic structure factor. For the tracer particle 
simulations we fix the initial position of a tracer particle 
at the center of the system with equally many particles 
to the left and right. The procedure above then has to 
be done separately for the particles to the left and right 
of the tracer particle. 

Step 7 and 8 constitute our rule for handling collisions 
between particles. These rules are based on the idea of 
having momentum conservation on average in collisions 
while at the same time fulfilling detailed balance. The to- 
tal momentum dissipated to the surrounding medium if 
particle n moves a distance I is l£ n . When the M+l parti- 
cles move the dissipated momentum is / Y^m^n Ca r " 
rying out the later move with a probability £„/ J^" £ m 
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makes the dissipated momentum independent, on aver- 
age, of whether collisions happen or not. Detailed bal- 
ance is satisfied, because the rate for the attempted move 
of the M + 1 particles is q n oc Thus the rate at 

which the particles actually move will be proportional to 

Vfn x WEl=!Cm = VELn^, and identical to 
the rate with which the reverse move (particle n ± M ini- 
tiating the reverse move of the M + 1 particles in the op- 
posite direction) occurs, since in equilibrium all allowed 
positions of the hard-core interacting particles are equally 
probable. Detailed balance is therefore satisfied by the 
above algorithm. 

In the simulations of HT systems we draw the friction 
constants £„ from the density [0 < a < 1] 

Q(Zn) = Mn 1 ^ for £ > £c (G5) 

and = for < £ c . The normalization con- 

stant is A = In order to generate a random fric- 

tion constant £ n according to the density above, in a 
standard fashion [34]], we set the cumulative distribution 
C(£n) = L n g{C)dS,' equal to a uniform random number 



r between and 1, C(£„) = r. Explicitly we then get 

Hn=Ul-r)- 1/a . (G6) 

In each simulation we determine a random friction con- 
stant £ n for each particle using Eq. (|G6[) . Then, using 
the fact that the free particle diffusion constant D n is 
related to the friction and hop rates as D n = q n a 2 /2 = 
ksT/^n, we get random hop rates: q n = 2k B T '/ '(£ n a 2 ) 
used in the simulations. We fix the normalization con- 
stant A by setting the average diffusion constant D = 
k B T j CV^n)^ = ak B T/[{\ + a)f c ] to unity in the 
simulations. 

In all simulations except those for the structure fac- 
tor S(Q, t) we have used a width a = 1 for the jump 
length distribution. For the structure factor simulations 
a = 0.25 was used due to a higher choice of density p. 
Ideally a should be as small as possible compared with 
the average distance between particles p _1 to approxi- 
mate the diffusive limit of moving particles in a spatial 
continuum. 
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